A Theoretical Search for the Optimum Giant Magnetoresistance 
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The maximum current-perpendicular-to-plane giant magnetoresistance is searched for in magnetic 
multilayers made of Co, Ni, and Cu with disorder levels similar to those found in room temperature 
experiments. The calculation is made possible by a highly optimized linear response code, which 
uses the impurity averaged Green's function technique and a 9-band per spin tight binding model. 
Using simulated annealing, hundreds of different configurations of the atomic layers are examined 
to find a maximum GMR of 450% in ultrathin Ni/Cu super lattices. 

PACS numbers: 75.70.Pa, 72.25.Mk, 72.15.Gd 
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The phenomena of the giant magnetoresistance (GMR) 
has been extensively investigated both theoretically and 
experimentally since its discovery, (fl, ^ This interest 
stems from both fundamental questions about magnetism 
and transport on the nanoscale and also from techno- 
logical applications like magnetic sensors. || From this 
body of work a basic understanding of the physics behind 
the GMR has emerged. When a magnetic field causes 
the magnetic domains in a magnetic multilayer or other 
GMR structure to change their relative orientation, the 
resistance changes because the bulk and interface scatter- 
ing rates depend on the electron spin orientation. While 
this is simple to explain, the actual details are quite com- 
plex. Experimentally, it is known that the GMR is very 
sensitive to atomic composition of the magnetic multilay- 
ers. In cases adding a single monolayer can change the 
magnetoresistance. W Film growth parameters can simi- 
larly change the GMR by an order of magnitude. H The 
sensitivity to atomic composition and film morphology as 
well as GMR samples consisting of many atomic layers 
has made the theory of the GMR challenging (for reviews 
see Refs. |andg). 

Nonetheless, much progress has been made in elucidat- 
ing the physics of the GMR in general and in explaining 
specific experiments. Thus, in this paper we take a dif- 
ferent approach. Based on our work and those of other 
groups, we assume that we have an accurate model of the 
GMR. H We then ask the question: what is the largest 
GMR one can obtain for a class of magnetic multilayers 
with experimentally reasonable disorder parameters? In 
other words, can one theoretically search the parameter 
space of a specific class of GMR structures and find the 
optimum or near optimum atomic configuration? This 
is a much more challenging problem than computing the 
GMR for a specific GMR structure because one must 
perform many individual GMR calculations in the search 
for the optimum magnetoresistance. Indeed there have 
been very few calculations like this in condensed mat- 
ter physics. Franceschetti and Zunger demonstrated the 
use of simulated annealing to optimize the band gap in 
semiconductor alloys and superlattices. By searching as 
much as I0 4 samples, they found the maximum band gap 



configuration to be very complex. Q Using a similar ap- 
proach, Ihiguez and Bellaiche optimized the ground state 
structural properties and electromachanical response of 
perovskite alloys. |fC| This article presents one of the first 
attempts to optimize a transport property. 

Optimizing the transport properties of nanostructures 
requires a delicate balance between accuracy and speed. 
To ensure that the optimization results are valid, the 
model has to include the band structure of the materials, 
the geometry of the nanostructure, and scattering effects 
caused by impurities, disorder, and temperature. On the 
other hand, the calculation of the transport property has 
to be fast enough that the optimization can be finished in 
a reasonable amount time. To solve these difficulties, we 
have developed a highly optimized parallel tight-binding 
code to compute electron transport in nanostructures us- 
ing the impurity averaged Green function technique and 
linear response theory. While the tight binding model we 
use is not as accurate as other techniques such as density 
functional theory, |Tf[ it does allow us to include real- 
istic band structure efficiently. Similarly, the impurity 
averaged Green's function technique allows us to include 
different kinds of interface and bulk scattering without 
worrying about the placement of individual defects. 

The tight binding Hamiltonian we consider is 



H = 



E 

r j;r' ,j' 



(1) 



where r and r' are sites on the lattice. The indices j 
and j' denote one of the one s, three p, or five d levels 
and also the spin orientation. Within a tight binding 
model, the disorder comes from variations in the matrix 
elements. In the impurity averaged green function tech- 
nique, these variations are treated statistically by aver- 
aging over an ensemble of different microscopic configu- 
rations of the disorder with the same macroscopic prop- 
erties, e.g., density of impurities. Denoting this aver- 
age by angular brackets, our unperturbed Hamiltonian is 
H Q = (H) . The deviation from this average for a particu- 
lar configuration of the disorder is given by 8H = H —H Q . 

Scattering is included through the choice of the self- 
energy. We use the self-consistent Born approximation 
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with uncorrelated disorder fluctuations. For energy ui the 
retarded and advanced self-energies, , are expressed 

in terms of the retarded and advanced green functions, 
G R / A as 

S^H = ((^r J;rJ ) 2 )G^M, (2) 

where the green function is in turn expressed in terms of 
the self-energy via 



[w - H Q - ?, R/a (lu)}G R/a (lu) = 1. 



(3) 



This last equation is a matrix equation, where the r and 
j indices have been suppressed. 

Within the linear response theory, the conductivity is 
computed via the Kubo formula. For the geometry with 
the current flowing perpendicular to the planes (CPP) it 
is necessary to include the vertex correction, T, which is 
given by the self-consistent equation, 



<(£tf rj;rj ) 2 }(G R (J + r)G 



where J is the current operator, 



-i(r - r')#, 



(4) 



(5) 



A conserving approximation for the conductivity at r is 



a a p{v) = jTi{J a (r)G R (J +T )G A 

- l -J a {v)G R J, 3 G R - ij Q (r)G A J /3 G A }, 



(6) 



where J(r) is the local current density operator at r. The 
sum of J(r) over position times the lattice unit cell vol- 
ume is equal to the current operator in Eq. (5). In the 
CPP geometry this conductivity will be constant because 
the current density is independent of position. 

Letting op be the conductivity when the magnetic lay- 
ers are aligned ferromagnetically and <Jj^ F be the conduc- 
tivity when adjacent magnetic layers are aligned antifcr- 
romagnetically, the giant magnetoresistance is defined as 
GMR = (op — °Af)/°Af- Note that in practice the 
magnetic layers will not necessarily be aligned antiferro- 
magnetically at zero field; however, our simulation does 
not compute the magnetic domain structure, but only 
the conductivity. In particular some of the ultrathin mag- 
netic layers which we find with large GMR could very well 
be ferromagnetically coupled. |l2| It is possible to achieve 
antiparallel alignment even with ferromagnetic coupling 
by going to a spin valve geometry with only two magnetic 
layers where one of the layers is exchange biased. p3{ 

An optimization study involves many GMR calcula- 
tions, therefore, the speed for each GMR calculation is 
crucial. The trace in Eq. (6) is over all positions and 
atomic levels; however, since the system is periodic, it 
can be broken up into a sum in k-space and a sum within 



a unit cell. Thus, the matrices are of order 9N for each 
spin, where N is the number of atoms in a unit cell. At 
first glance, the matrix products and inverses in Eqs. (3), 
(4) and (6) are all 0(N 3 ) operations, where N is the num- 
ber of atoms in the unit cell. In addition, both Eqs. (3) 
and (4) are self-consistent equations. Performing opti- 
mization with a self-consistent 0(N 3 ) implementation is 
estimated to take about 5 CPU-years, which is impracti- 
cal. To overcome this, we have developed an O(N) algo- 
rithm to solve Eqs. (3) and (6), and an 0{N 2 ) algorithm 
to solve Eq. (4) directly without the need for iteration. 
For N = 50, this algorithm is more than 100 times faster 
than the 0(N 3 ) algorithm. Running on a Linux cluster 
of five dual 866MHz Pentium III CPU workstations, it 
takes about 20 minutes for each GMR calculation and 
about two weeks for an optimization. 

The optimization technique we use is simulated an- 
nealing. Since not every configuration qualifies as a 
GMR structure, the optimization is performed in the 
constrained subspace consisting of a unit cell with two 
magnetic layers separated by two non-magnetic spacer 
layers. After an arbitrarily chosen initial unit cell, each 
subsequent configuration is generated by one of the fol- 
lowing Monte Carlo moves: (i) inserting a monolayer, (ii) 
removing a monolayer, or (iii) changing the composition 
of a monolayer. If the GMR of the new configuration 
is higher than the previous configuration, then it is ac- 
cepted. On the other hand, if the new GMR is lower, 
it is only accepted with a probability of exp(AGMR/T), 
where AGMR is the change in the GMR and T is the 
simulated annealing temperature. The process is contin- 
ued as the annealing temperature is gradually lowered 
until there is no further change in the configuration and 
a final near-global maximum in the GMR is attained. 

In the following, we show results of a study that looks 
for the optimal configuration of superlattices made with 
Co, Ni, and Cu in conditions similar to GMR experi- 
ments at room temperature. To make the study more 
manageable, we only consider fee (111) lattices. The 
minimum thickness of each layer is limited to two mono- 
layers in order to avoid complications such as pin holes. 
The tight binding parameters in Eq. (1) are obtained 
from fits to density-functional calculations where up to 
second-nearest-neighbor hoping energies are included. 

Two types of scattering are included: spin-independent 
scattering and spin-dependent scattering. The spin- 
independent scattering represents the structural disorder 
and is included in the calculation by uniformly shifting 
all the orbital levels at a site: {(5H r j- r j) 2 ) — 0.2eU 2 , 
which corresponds to a Cu resistivity of 3.1 /iilcm, about 
twice the value of clean Cu at room temperature. The 
spin-dependent scattering is assumed to be due to inter- 
diffusion of atoms near the interface. It is modeled by 
the shifting of the (i-levels only, since the d-levels are the 
primary difference in the tight binding parameters for 
these elements. The shift in the d-levels are determined 
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FIG. 1: GMR as a function of the simulated annealing step. 
The four marked configurations and their GMR values are: 
(a) [C03/CU3/C03/CU3] 36%, (b) [Ni 2 /Cu 3 /Ni 2 /Cu 2 ] 221%, 
(c) [Nia/Cua/Nia/Coi/Cua] 38%, and (d) [Ni 2 /Cu 2 /Ni 2 /Cu 2 ] 
450%. Notice that the GMR is very sensitive to the configura- 
tion. Both (b) and (c) differ from (d) by only one monolayer. 



by fits to 4-atom supercell density-functional calculations 
of an impurity in three host atoms. We then assume an 
interface mixing percentage of 10% in each of the two 
monolayers on either side of the interface. The value of 
{(SH r ,j- r ,j) 2 ) for the majority and minority rf-states is 
equal to 0.29 and 0.84 eV 2 at the Co/Cu interface, and 
0.03 and 0.11 eV 2 at the Ni/Cu interface. Using these pa- 
rameters, the GMR for a [C03/CU3/C03/CU3] super lat- 
tice is 36%, which is in the range of observed CPP-GMR 
at room temperature. Jl4| Our spin-dependent scattering 
does not include spin flip scattering, which can be im- 
portant for short spin flip scattering lengths. Jl5| 

The giant magnetoresistance for a simulated annealing 
run is shown in Fig. 1. The initial configuration has been 
chosen to be the superlattice [C03/CU3/C03/CU3], with 
a unit cell made of four layers, each consisting of three 
monolayers of either Co or Cu. The particular choice of 
initial configuration does not effect the final result for this 
calculation. As the annealing temperature T decreases 
from its initial value of 60%, the GMR becomes larger. 
The GMR does not increase monotonically because a su- 
perlattice that reduces the GMR is accepted with a prob- 
ability set by the Boltzmann factor, which is crucial for 
escaping from local maxima in the GMR. As the anneal- 
ing temperature is lowered, the chance of leaving a local 
maximum becomes smaller and smaller. At step 89 in 
this run, the GMR jumps to 450% and no other config- 
urations are accepted. The atomic configuration for the 
highest GMR (d) is [NIa/Ctfe/Nia/Cifc]. This is a very 
surprising result because most room temperature studies 
to date have GMR values less than 100%. IM 111 111 From 



FIG. 2: Histogram of the GMR in the configuration space for 
600 randomly generated configurations. The average GMR is 
22%, and the standard deviation is 36%. About 97% of the 
configurations are located near the main peak with GMR val- 
ues less than 100%, while 3% of the configurations are many 
standard deviations away from the main peak with GMR val- 
ues greater than 100%. These large GMR configurations are 
mostly ultrathin Ni/Cu superlattices. 



Fig. 1 one can see that the GMR is very sensitive to the 
atomic configuration. For example, superlattices (b) and 
(c) have GMR values of 221% and 38%, respectively, yet 
they differ by only one monolayer from the optimal GMR 
of superlattice (d), 450%. 

To get a better picture of the GMR in configuration 
space, we have calculated the GMR for 600 randomly 
generated configurations. The result is shown in Fig. 
2 as a histogram. The average GMR is 22%, and the 
standard deviation is 36%. There are two groups of 
configurations. About 97% of the configurations are lo- 
cated near the main peak. They have a GMR less than 
100%. On the other hand, about 3% of the configura- 
tions have GMR many standard deviations away from 
the main peak. They scatter in the range from 100% to 
450%. Most of the large GMR configurations are ultra- 
thin Ni/Cu superlattices. Only two of the large GMR 
configurations contain Co. Since the large GMR config- 
urations are rare, they may have been overlooked. 

What causes the GMR to be large in some multilayers 
but not others? The dominant source of scattering in 
these transition metals is due to s-d hybridization. The 
rf-levels in the majority spin of Ni, Co, and Cu are below 
the Fermi level; however, there is still a finite d-density 
of states at the Fermi level due to both s-d hybridization 
and disorder. In particular disorder can substantially 
smear the rf-density of states, increasing the value at the 
Fermi level. || Since the rf-levels in Ni and Cu are more 
well matched than those in Co and Cu, the disorder due 
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In this paper we have performed a simulated annealing 
search for the optimum or maximum giant magnetoresis- 
tance by varying the atomic composition of multilayers 
composed of Ni, Co, and Cu. The disorder parameters 
in our calculation were chosen to simulate room temper- 
ature experimental conditions. After hundreds of GMR 
calculations for different atomic configurations, a surpris- 
ingly large GMR was found for ultrathin Ni/Cu multi- 
layers. The origin of this large GMR of 450% was found 
to be the fact that the (i-levels in Ni and Cu are relatively 
close in energy. Thus, although we have made a number 
of approximations in order to perform our calculation in 
a reasonable amount of time, the large GMR appears to 
be a robust effect. We believe that this work is part of a 
large class of problems in tailoring transport properties 
of nanostructures to optimize an observable such as the 
magnetoresistance. 
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sions. This work is supported by DOD/AFOSR F49620- 
1-0026, NSF DMR9357474, and the Research Corp. 



FIG. 3: Comparison of [Ni 2 /Cu2/Ni 2 /Cu2] (circles) and 
[C02/CU2/C02/C112] (crosses) superlattices. (a) At the Fermi 
level the cf-density of states for the majority spin when the 
magnetic layer aligns ferromagnetically, rf-DOSp aj , is higher 
in the Co/Cu sample, leading to higher scattering, (b) The 
conductivity for the same spin channel, o"JJ laj , in Co/Cu is 
thus lower than that in Ni/Cu. Since near the Fermi energy 
the conductivities for the minority spin and the antiferromag- 
netically aligned magnetic layers are relatively small and flat 
compared with a™ ai , the GMR (c) resembles a^ a3 . 

to interface mixing is smaller in Ni/Cu multilayer than 
in a Co/Cu multilayer. Thus, the GMR is greater in a 
Ni/Cu multilayer for the same level of interface mixing 
because scattering via the s-d hybridization is less. 

To see this more clearly, in Fig. 3 we take a represen- 
tative Ni/Cu sample and a representative Co/Cu sam- 
ple. The two samples have the same bulk disorder, but 
the Ni based multilayer actually has more interface mix- 
ing (10%) than the one containing Co (5%) to allow the 
curves to fit clearly on the same plots. As seen in Fig. 
3(a), the d-density of states is smaller in the Ni sample. 
The d-density of states correlates well with the larger con- 
ductivity shown in Fig. 3(b) and the larger GMR shown 
in Fig. 3(c). The fact that there is a simple explanation 
for this large GMR indicates that it is a robust effect and 
not strongly dependent on the many assumptions which 
we have made. Also, it says that for comparable coupling 
between layers and disorder, the Ni/Cu multilayers will 
have a larger GMR than the Co/Cu ones. 
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